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Abstract In this article we investigate the permeability of a porous medium as 
given in Darcy’s law. The permeability is described by an effective hydraulic pore 
radius in the porous medium, the fluctuation in local hydraulic pore radii, the 
length of streamlines, and the fractional volume conducting flow. The effective 
hydraulic pore radius is related to a characteristic hydraulic length, the fluctuation 
in local hydraulic radii is related to a constriction factor, the length of streamlines 
is characterized by a tortuosity, and the fractional volume conducting flow from 
inlet to outlet is described by an effective porosity. The characteristic length, 
the constriction factor, the tortuosity and the effective porosity are thus intrinsic 
descriptors of the pore structure relative to direction. We show that the combined 
effect of our pore structure description fully describes the permeability of a porous 
medium. The theory is applied to idealized porous media, where it reproduces 
Darcy’s law for fluid flow derived from the Hagen-Poiseuille equation. We also 
apply this theory to full network models of Fontainebleau sandstone, where we 
show how the pore structure and permeability correlate with porosity for such 
natural porous media. This work establishes how the permeability can be related 
to porosity, in the sense of Kozeny-Carman, through fundamental and well-defined 
pore structure parameters: characteristic length, constriction, and tortuosity. 


1 Introduction 


Understanding the process of flow through porous media is of great importance 
in many fields, including petroleum engineering and hydrology. Slow viscous fluid 
flow in porous media is traditionally described by Darcy’s law, a proportional 
relation that links the fluid discharge Q to an applied piezometric head (hydraulic 
head) difference Ah: 


Qjj.As 

ApgAh^ 
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where A is the cross-sectional area of the porous medium, As is the length of 
the porous medium in direction of applied head difference, n is the constant fluid 
viscosity, p is the constant fluid density, g is acceleration due to gravity, and k 
is a constant for the porous medium called (intrinsic) permeability |13lF7l[T^. For 
simplicity, pgh will also be called the piezometric head throughout this article. 

A fundamental question in flow through porous media is how the permeability 
k can be related to the porosity (j) through well-dehned parameters of the pore 
structure. A well-known relation of this kind is the semi-empirical Kozeny-Carman 
equation. It was first proposed by Kozeny |24| as 


k cx T 




( 2 ) 


where dw is an effective grain size, and r is a tortuosity of the porous medium 
describing the relative difference between the microscopic (interstitial) head gra¬ 
dient along the streamline and the macroscopic head gradient. Kozeny derived 
his equation assuming that the porous medium could be viewed as a bundle of 
streamtubes |24| . 

Carman noted that linking the microscopic fluid velocity to the Darcy velocity 
for the porous medium involves scaling with the factor r m- Carman therefore 
modified Kozeny’s Eq. (ED by multiplying with the tortuosity r m- 

OC (3) 


For a monodisperse sphere pack we have dw = 6(1 — 4>)/So, where 5'o is the spe¬ 
cific surface area EH- This leads to a more general form of the Kozeny-Carman 
equation: 

k = cor^r^0, (4) 

here r/j = (j)/So is the (mean) hydraulic radius and Co is a coefficient called Kozeny’s 
constant miMm- 

The hydraulic radius r/j is assumed to represent an effective pore radius of the 
porous medium. It is a purely geometric length which does not take into account 
the effect on permeability from pore size variation or connectivity. Others have 
proposed length-scales more suitable for permeability description, among these 
the smallest pore along the most conductive percolating pathway HD, nuclear 
magnetic resonance relaxation time [6], or grain size distribution m- Johnson 
et al. |2n| suggested a length by weighting with the electric field, thereby related 
to (electrical) transport and thus a dynamical length in contrast to the geometri¬ 
cal hydraulic radius r^. This length has been shown to be a better permeability 
descriptor than the hydraulic radius |36| . however, there is no fixed relation 
between electrical conductance and fluid flow E2- A dynamical length linked to 
fluid flow instead of electrical conductance was introduced in [8]. This length was 
derived from the microscopic hydraulic conductance [8117] . thus descriptive of fluid 
flow in porous media. 

For a porous medium, the (hydraulic) tortuosity is a measure of the microscopic 
flows deviation from the direction of the applied piezometric head difference, re¬ 
flected by the length of the microscopic streamlines Email]. Tortuosity has been 
defined as As/se, where As is the length of the porous medium in the direction 
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of applied piezometric head difference, and Se is the effective streamline length [71 

[22] . There is ambiguity associated with the derivation of Se [71(2^ ■ however most 
formulations are a weighted average of the streamline lengths m- 

Furthermore, the constricting and expanding nature of pore channels converges 
and diverges the streamlines, which leads to variation in fluid velocity along the 
streamlines and decrease the permeability. This effect has been treated for simpli¬ 
fied porous media m, for more complex material the standard deviation of the 
cross-sectional area has been used as an estimate |3l] , while for electrical conduc¬ 
tance the effect of pore channel variation has been described in general j^. We 
follow Ref. |15| and use the term constriction factor to account for this porous 
medium property. Note that the term constrictibility is common when considering 
effective diffusion in porous media, see e.g. Ref. isil. The effect on transport from 
constrictions is frequently lumped together with the effect from tortuosity, see e.g. 
Refs. [SBlIHlIT] . However, the pore structure effect on permeability from these two 
distinct geometrical properties can be separated for any porous medium, as shown 
in this article. 

The (geometrical) porosity (p = QjV is the fraction of pore space f? in the 
porous medium of total volume V. The pore space Q is sometimes interpreted as 
the connected pore space, thus both the permeability and porosity vanish at a 
percolation threshold at a finite geometrical porosity m- Other authors have also 
excluded dead-end pores [TlfTS] . or considered an effective porosity based on the 
streamlines connected both to the inlet and the outlet of the porous medium [H 

[23] . 

A large body of literature exists on the relation between macroscopic trans¬ 
port properties and pore structure; the interested reader is referred to m and 
[7] for reviews of numerous relations. The recent progress in computing and imag¬ 
ing provides tools for advances on the topic. A number of methods for statisti¬ 
cally EllMlIISl or process based m reconstruction of three-dimensional porous 
media from two-dimensional thin section images has been proposed. Advances 
in micro-computed tomography (mCT) have made it possible to directly image 
many types of natural porous media with sufficient resolution to represent the 
three-dimensional pore structure [35li|. For flow and transport simulations net¬ 
work analogs have been widely used to represent the pore space [I7|[II1[29|, while 
the full microscopic pressure and velocity fields can also be computed for both 
idealized models [251IM1I^ and for natural porous media such as reservoir rocks 

mEi- 

The aim of this work is to derive a comprehensive relation for porous media be¬ 
tween the permeability and porosity from detailed pore structure information. In 
contrast with existing relations the permeability will be fully defined by separable 
descriptors of the pore structure without introducing free parameters or constants. 
In Sect. [2| we introduce a microscopic permeability factor which describes the lo¬ 
cal contribution to the effectiveness of the pore space to conduct fluid flow. This 
permeability factor is decomposed into streamlines of the fluid flow, and further 
factorized into distinct contributions from characteristic length, constriction and 
tortuosity in Sect. [3| In Sect. Hj we integrate the characteristic length, constric¬ 
tion and tortuosity from individual streamlines into effective parameters, all pore 
structure descriptors. In Sect. [5| we use the Hagen-Poiseuille equation to demon¬ 
strate our approach on idealized porous media, while the methodology is applied 
to network analogs of Fontainebleau sandstone data in Sect. [6| 
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2 Microscopic (Interstitial) Permeability 

Consider a porous medium V of length As in direction of an applied piezometric 
head difference Ah, consisting of matrix and pore space f2 C V filled with an 
incompressible fluid. At the microscopic (interstitial) scale, a slow (creeping) flow 
is governed by the Stokes equation supplemented by the continuity equation: 


/iV^u = pgVh, 


(5a) 

(5b) 


V • u = 0, 


where u is the microscopic fluid velocity, and h is the microscopic piezometric head 
msM- In the following we will refer to Eqs. © simply as the Stokes equations. 
Throughout this article the fluid flow is assumed to be governed by the Stokes 
equations. There are no time derivative terms in the Stokes equations; hence a 
constant piezometric head difference Ah implies a steady-state flow. 

Let S denote the set of all streamlines 5 connected both to the inlet and the 
outlet of the porous medium, and let f7s = {x € 5 € S} be the subset of f? where 
the fluid flows from inlet to outlet. The effective porosity is (j)s = SlsjV Eina- 
Due to the linearity of the Stokes equations, the streamlines 5 are independent of 
the magnitude of applied piezometric head drop Ah and the constants p and p, 
thus Ds and (jjg are only dependent on pore structure and direction of the applied 
piezometric head drop. Note that for dead end pores we might have pore space 
that is not in still the fluid velocity might be non-zero, called reentrant flow 
in Ref. [H]. 

In our system the piezometric head difference —pgAh is the potential that 
drives the fluid through the porous medium, hence the rate of applied energy is 
—pgAhQ. The potential driving the microscopic fluid flow is the piezometric head 
pgh, and the rate of work done by this piezometric head potential is given by 
—pg\7h ■ u. 

Applying the divergence theorem, and invoking that the fluid velocity u is a 
solenoidal vector field from the continuity equation Eq. (I5bl) . we obtain 



pg{houtQ - hinQ) = pgAhQ, 


( 6 ) 


where n is the outward pointing unit normal field of the boundary SOs, hin is the 
piezometric head at inlet, and hout the head at outlet. Here we use that u • n = 0 
except at inlet and outlet. 

The permeability as given by Darcy’s law in Eq. © can then be expressed as 
follows: 



We will denote the integrand in Eq. © as the microscopic permeability factor'. 



(8) 
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The microscopic permeability factor k is then the rate of work done by the piezo¬ 
metric head potential —pgVh ■ u multiplied by the factor pAs^/{pgAh)^, however 
it is a constant only dependent on the pore space n C V and the direction of the 
applied head difference Ah, in contrast to —pgVh ■ u. 

We derive the effective permeability factor by integrating k over Qs '■ 

Since the microscopic permeability factors k are constants only dependent on the 
pore space and the direction of the applied head difference Ah, so is the macro¬ 
scopic permeability factor Ks. Moreover 


k = Ks4)s- (10) 

We have thereby divided the permeability into two factors; the effective porosity 
(fs yielding the pore space fraction where fluid flow from inlet to outlet, and the 
permeability factor k.s yielding the effectiveness of the pore space fls to conduct 
fluid flow. Both (fs and Ks are dependent on direction, leading to the anisotropy 
of the permeability. 

Following the same arguments as above we also have 


k = kef. 


( 11 ) 


where 


This implies that 


K = 





ndV = 0, 


( 12 ) 

(13) 


even though the permeability factor k might be non-zero in 17\17s. Hence reentrant 
flow does not contribute to the effective permeability factor. In the following we 
work with streamlines connected to inlet and outlet, thereby excluding the part of 
the pore space containing reentrant flow in addition to stagnant parts. 


3 Streamline Decomposition 

In this section, we show how the permeability factor can be decomposed onto 
streamlines, and furthermore how the permeability factor for each streamline can 
be segmented into parts describing the hydraulic conductance, the constrictions 
and the tortuosity along this streamline. 

We can discretize the space 17s into a disjoint union of simply connected spaces, 
such that each streamline is fully contained within one simply connected space. 
Using the continuity equation given by Eq. there exist scalar functions A 

and X such that VH x VX = u [71131 . The scalar functions A and X represent two 
families of stream surfaces whose intersections are the streamlines. 

Every point in a; £ 17s is uniquely described by the streamline 5 passing through 
X, and the distance s along S from inlet to point x. The streamline S = 5(71 = 
A, X = x) is the intersection of the surfaces given by the constants A and x- This 
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is similar to a Lagrangian frame of reference, however we use distance instead of 
time to distinguish points on a streamline. A change of variables from the usual 
Cartesian coordinates {x, y, z) to the streamline coordinates (A, X, s) gives the 
Jacobian 


S{A,X,s) 

S{x,y,z) 


(VA X VA) . Vs = u - - = u. 

u 


(14) 


By vector calculus identities we have u = VA x VA = V x AVA. Using Stokes’ 
theorem, the fluid discharge Q through a streamtube bounded by four stream 
surfaces represented by the constants Ai,A 2 ,Xi and X ‘2 is then 


Q= f f u ndS= f AVA.lds = (A 2 -Ai)(x 2 -Xi), (15) 

J Js Jss 

where 5 is a cross section of the streamtube, and 1 is the unit tangent to the surface 
boundary 3S. In the last equality we use that either A or A is constant for each 
of the four line segments in SS. 

We can now define a permeability factor for the individual streamlines. Starting 
with Eq. we have 


Ks = 


t.L 

tJIL 

T.l 


KdV = 


Qt 


JJL 


-dsdXdA 


SiA,X) 


Vh ■ u f As \ ‘ 

-m -I —^ 1 


5(yl.X) 


-ypgAh 


As y 

pgAh) 


dQs = 


\pgAhJ 
1 


Us 


dsdXdA 


K{S)dQs, 


(16) 


where k{S) = —pAs^/{pgAh) is the permeability factor for the streamline S. For 
the fourth equality we use that V/i-u/u = Sh/5s, and that the total head difference 
along a streamline is equal the applied head difference Ah. The inhnitesimal fluid 
discharge for the infinitesimal streamtube given by dX and dA is denoted by dQg, 
where dXdA = dQg from Eq. (US. 

Note that k(S) is a constant, and that 

= W ! <S)dQs = 7 ^ k ( 5 ) / dQs = -^^( 5 ). ( 17 ) 

Rewriting the expression for k( 5), we obtain: 

k(5) = f - ^ — dsf^Ah [ ^ ds) , (18) 

pgAh \ls J Js P^Vh-u jgVhu J 


where Is is the length of the streamline 5. In the following, we will link the 
permeability factor for a streamline to descriptors of the pore structure: namely, 
tortuosity, constriction and hydraulic conductance; all represented in Eq. m- 
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3.1 Tortuosity 

The tortuosity of the streamline 5 is given by r(5) = As/lg, i.e. the length of the 
porous medium divided by the length of the streamline 0122]. Due to the linearity 
of the Stokes equations, Eqs. m, the streamline 5 is independent of the magnitude 
of applied piezometric head drop Ah and the constants p and p, thus t{S) is only 
dependent on pore structure and direction of the applied piezometric head drop. 

For smaller tortuosity r(5) the fluid needs to travel longer distance, and more 
applied head potential is expended due to transport distance. This increase in 
energy expenditure is reflected in the smaller factor t(S)^ = (As/lg)^ in Eq. (I18II . 
Longer travel distance for the fluid decreases the effectiveness of the pore space to 
conduct flow. 

3.2 Constriction Factor 

For a straight circular pore channel of length L with cross-sectional area A{x) at 
point X, the degree of variation in cross-sectional area can be measured by the 
constriction factor 



(19) 


( 20 ) 


corresponding to dehnitions introduced in Refs. [T5l[9|. For the second equality we 
assume the fluid flow is described by the Hagen-Poiseuille equation (see Eq. (1331) 1. 
thus Q/{pgVh{x)) oc A{x)^. When the fluid is incompressible, the total discharge 
Q must be constant through all pore channel cross-sections due to mass-balance, 
yielding Eq. (1201) . For porous media in general, the cross-sectional area A(x) used 
in Eq. (ED is not straight-forward dehned. As seen in Sect. 15.41 C represents the 
reduction in permeability due to the variation in cross-sectional area. 

Following Ref. [9| we propose a (hydraulic) constriction factor for streamline S 
by replacing the head gradient V/i in Eq. (1201) with the head derivative 5h/5s = 
V/i • u/u along the streamline: 



( 21 ) 


As with the tortuosity r(5), the constriction factor C{S) is only dependent on 
pore structure and direction. 

When the fluid flows through a constriction, the head derivative 5h/Ss along 
the streamline increases. A large variation in pore size along the streamline then 
translates into a large variation in the head derivative. The constriction factor C{S) 
thus relates to the constricting and expanding nature of the pore space along the 
streamline 5, or equivalently the converging-diverging set of streamlines around 
streamline 5. For a larger constriction factor C(5) the effectiveness of the pore 
space to conduct flow is reduced. 
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3.3 Hydraulic Conductance 


The microscopic hydraulic conductance is given by 


B = - 


pgVh ■ u’ 


( 22 ) 


and is related to the pore size and shape, and the location in the pore [81171. 
Following Eq. m, the hydraulic conductance for a streamline S is 


B(S) = [ B-ds = [ - ^ — ds. (23) 

Js ^ Js pg^h-u ^ 

Observe that B{S) represents the second term of Eq. (1181) . Also note that B{S) 
is dependent on both the magnitude of the applied piezometric head drop pgAh 
and the viscosity p, in addition to pore structure and direction of the applied 
head drop. This is in contrast with the hydraulic conductance B, tortuosity r(5) 
and constriction factor C(5), which are only dependent on pore structure and 
direction. 

With the formulations above, Eq. m can now be rewritten as follows: 


k{S) 


C(S) 


(24) 


The permeability factor for an individual streamline is then expressed by descrip¬ 
tors of the pore structure. 


4 Effective Permeability 

In this section we show how the permeability factor ks can be segmented into a 
characteristic length, a constriction factor and a tortuosity by averaging over the 
streamline values. These are strictly related to fluid flow, based on the solution of 
the Stokes equation inside the porous medium. 

The effective hydraulic conductance is found as the volume-weighted average 
of the hydraulic conductance B [81171 : 

Bs = ^ f BdV. (25) 

J f2s 

Since B is only dependent on pore space and direction, so is Bs. From 

= 7T / B{S)dQs, (26) 

JOs JS 

we have a correspondence between the volume integral of B and the streamline 
integral of B{S). 

We define the characteristic (hydraulic) length as L/j = y/8Bs to represent the 
effective hydraulic pore radius of the porous medium. Note that the characteristic 
length scales linearly with the size of the porous medium, as desired for a char¬ 
acteristic length. As seen in Sect. 15.21 for a porous medium consisting of parallel 
circular tubes of radius r and length As, where these tubes connect the opposite 
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sides of a cube of side length As, then = r, and k = 4>si^s = (psBs = 4>sL\/9> as 
desired for such a medium m- 

Consider a porous medium for which C{S) = 1 for all streamlines S, e.g. a 
single tube of constant cross-sectional area. Following Eq. dl. for such a porous 
medium we desire for a tortuosity Ts to give k = (j>sBs Tg |15lF7|. Then ks = BsTs 
by Eq. (Uni), thus invoking Eq. (El gives 

= r RtLn [ (27) 

jsB{S)dQs Js 

Hence the tortuosity squared t] is a weighted average of the streamline tortuosity 
squared. Note that the tortuosity is only dependent on the pore space Qs 
and the direction of applied piezometric head drop, it is dimensionless and scale 
invariant. 

The tortuosity is commonly formulated as a weighted average of the streamline 
lengths [Til[23] , therefore let the tortuosity = 1/ J^w{S)dQs J^T°‘{S)w{S)dQs 
be another weighted average of streamline tortuosity r(5), now to a power a. Con¬ 
sider a porous medium for which C{S) = 1 and r(5) is constant for all streamlines 
S. We still want k = (j)sBsT°‘, thus f“ = , which yields 

^ f = 7“ = = T^{S), (28) 

Js‘w(S)dQs Js 

therefore a = 2. If we desire a streamline decomposition to hold for Bgf^ = 
J^K(S)dQs, then {J^B{S)dQs/ J^w(S)dQs)w{S) = B{S), yielding unique of 
the form described by r“. 

Using Eqs. (fTTI) . (El, (1^51) and (1771) . we have 

Bst! = 7 ^ B{Sy{S)dQs = 7 ^ ^ K{S)C{S)dQs = C{S)dQs. (29) 

By factoring out the hydraulic conductance Bg and tortuosity , the remaining 
contribution to the permeability factor Kg is 

Cg = ^ f C{S)dQs, (30) 

^ Js 


where Cg is denoted the (hydraulic) constriction factor. The constriction factor is 
also only dependent on the pore space Us and direction, it is dimensionless and 
scale invariant. While the hydraulic conductance represents an effective hydraulic 
pore radius, the constriction factor represents the fluctuation in hydraulic pore 
radii. 

From Eq. (1291) and Eq. (1301) we have 


_ Bgri _ tUI 

Cg SC's ■ 

Combining Eq. (HOI) with Eq. (1311) then gives: 

k = KgJ>g = ll^ = ^l^ 

Os oOs 


(31) 


( 32 ) 
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We thus have a full description of the porous medium permeability by pore struc¬ 
ture related parameters. 

The permeability factor ks gives the effectiveness of the pore space 17s to 
conduct flow. This effectiveness is reduced by longer flow paths given with a smaller 
tortuosity ts, more variation in pore size along the flow paths described by a larger 
constriction factor Cs, and smaller pores reflected by a smaller characteristic length 
Lfj^. Note that these factors are dependent on direction in addition to the pore 
structure, which leads to the anisotropy of the permeability. 


5 Single Tube Example 


Capillary bundle of tube models have a wide use as simplified representations 
of porous media. The single tube examples in this section illustrate such model 
representations. 

Consider a straight cylindrical tube with constant cross-section of radius R. If 
the length of the tube is much larger than the radius, then the flow inside the tube 
is approximated by the Hagen-Poiseuille equation 


Q 


■wR'^pgVh 

8)1 


Moreover, the flow velocity is given by 


(33) 


u{r) 


{R^ -r'^)pgVh 

4m 


(34) 


where r is the distance from the center of the tube [T]. In the equation above and 
subsequently in this section, V/i is also used to denote the scalar value —||V/i||. 
In the following subsections, the Hagen-Poiseuille equation is used to demonstrate 
our theory introduced above. 


5.1 Permeability Factor 


Consider a straight cylindrical tube inside a cube of side-length L, and with an 
applied piezometric head difference Ah over two opposite sides of the cube. Let 
the tube be of length L and aligned with the applied head, then the piezometric 
head gradient inside the tube is Vh = AhjL. 

Combining Darcy’s law with the Hagen-Poiseuille equation, Eqs. m and (|33ll. 
the permeability of the cube is given by: 


QpL _ nR'^ 
UpgAh ^ SU' 


(35) 


Using the flow velocity a distance r from the tube center as given by Eq. pill, 
and the piezometric head gradient Vh = Ah/L, the permeability factor described 
in Eq. (jS]) is «:(r) = {R^ — Taking the volume-weighted average then gives 

the effective permeability factor 


= I (36) 

This gives = ■kR'^/{8L^), which is equal to the result from Darcy’s law in Eq. 
(1351) . hence our results are consistent with Eq. (nni. 
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Fig. 1 Idealized porous medium with tor- Fig. 2 Idealized porous medium with a 
tuosity different from 1. single constriction. 


5.2 Characteristic Length 


Let our porous medium and applied piezometric head be as above. The perme¬ 
ability can be calculated using the individual contributions from the characteristic 
length, constriction, and tortuosity, as given by Eq. (I32II . The tortuosity and con¬ 
striction are equal to 1 when the fluid flow inside the tube is described by the 
Hagen-Poiseuille equation, while the hydraulic conductance is B{r) = {R^ — r^)/4. 
The volume-weighted average of B(r) is then equal to the integration of «:(r) in Eq. 
(1361) . therefore Bg = /S. This gives a characteristic length = ^/8Bs = R. The 

characteristic length thus equals the radius of the tube, as desired. Since ts = 1 
and Cs = 1, we have 


8Cs 8L2 ’ 


(37) 


which agrees with the result from Darcy’s law in Eq. (ESI, hence our results are 
consistent with Eq. (13211 . 


5.3 Tortuosity 


Let again our porous medium and applied piezometric head be as above, except 
that the length of the straight circular tube with constant cross-sectional area is 
s > L, as shown in Fig. [T] Then the piezometric head gradient inside the tube is 
Vh = Ahjs. 

Combining Eqs. © and (|33| gives the permeability as: 


QpL _ 7ri?^ 
L'^pgAh 8Ls 


(38) 


Note that a larger s gives a smaller permeability k, while s = L gives a permeability 
equal to Eq. ([35|. 

The calculations for the characteristic length in Sect. I5.2l still hold, so L/j = R. 
The constriction factor Cs is equal to 1, while each streamline has length s, thus 
Ts = t{S) = L/s. We then have 


8Cs 8Ls ’ 


(39) 


which is equal to the result from Darcy’s law in Eq. (13811 . and consistent with Eq. 
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Note that by changing the orientation of the tube shown in Fig.[T]the tortuosity 
and the permeability will change, while in this idealized case the characteristic 
length and constriction factor stay constant. In this case the anisotropy of the 
permeability is captured by the tortuosity. 


5.4 Constriction 


We will now use an idealized porous medium to investigate the effect of con¬ 
striction. Consider a porous medium consisting of two tube segments in sequence 
aligned with the applied piezometric head, as depicted in Fig. [2] The two tube 
segments both have length I ,/2 inside a cube of side-length L, while the radii of 
the two segments are Ri and i? 2 - When assuming Ri<^ L12 for i = 1,2, the flow 
inside the tube segments can be approximated by the Hagen-Poiseuille equations 
and the tortuosity can be approximated as ts = 1 . 

Invoking Eq. (I33II we can show that the discharge is Q = 
— ( 7 r 7 ?fi? 2 P(?V/i)/( 4 p(i?f + Rt))- Using Darcy’s law, Eq. ([1]) , we then obtain 


Q^L _ nRfR^ 
L'^pgAh ~ AL^{Rf + R^)' 


(40) 


From Eq. (I33II we have V/ij = —8fj,Q/{TvRfpg), and from Eq. (12111 we derive the 
constriction factor as 


Cs = 



HRf + Rir 

4 RfR^ 


(41) 


For each tube section we have = Rf/S, which gives Bs = (i?f-|-i?|)/( 8 (i?f-|-i? 2 )). 
Since Tg = 1, we have: 


^ TvRfR^ 

Cs 4L2(Bf-hB|)’ 


(42) 


consistent with Eq. (liUll . Note that a larger difference between Ri and i ?2 gives a 
larger constriction factor Cs , which then implies a lower permeability, while Cs = 1 
when Ri = R 2 , as desired. 

We will now revisit the more general constriction example from Sect. 13.21 where 
we considered a tube with cross-sectional area A{x) for x £ [0,L]. We still as¬ 
sume that the flow inside the tube is approximated by the Hagen-Poiseuille equa¬ 
tions and that the tortuosity can be approximated as ts = 1. Then Vh(x) = 
—8pQTv/{A{x)^pg), from Eq. ([1]) we then have 


k= --. (43) 

The constriction factor 

equals Eq. (US in Sect. 13.21 For the cross-section at point x we have B{x) = 
ri(a;)/( 87 r), thus Bs = f A(x)^dx/(8Tr f A(x)dx). 
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Fig. 3 Visualization of the network representation of the pore space of Fontainebleau sand¬ 
stone mCTc in Tabel fTI Balls represent pore bodies, and sticks represent pore throats. This 
network has 7413 pore bodies and 14254 pore throats. 


Now assume another porous medium with equal pore volume and with a con¬ 
stant cross-sectional area. Then the cross-sectional area is A = (l/T) f A(x)dx, 
Bs = A/(9,-x) and k = A^/{L^8n). When factoring out the hydraulic conductance, 
the reduction in permeability due to the varying cross-sectional area is 


k/Bs ^ J_ 

k/Bs Cs ’ 


(44) 


hence the reduction equals the inverse of the constriction factor. 


6 Fontainebleau Rock Example 

We next turn to natural porous media, such as given by micro-CT (microtomog¬ 
raphy) images and rock models of Fontainebleau sandstone. Using the e-Core soft¬ 
ware m we generated three-dimensional rock models of Fontainebleau sandstone 
with porosities ranging from 8 to 26%. We used the exact same grain packing for 
all models, while we changed the amount of quartz cementation to achieve the 
variation in porosities. The rock modeling process is described in detail in Refs. 
[30119]. The micro-CT images and model sample size were 2.7 mm cubed with a 
resolution of 5.7 pm. 

We extracted network analogs using the e-Core software [T^. The software 
extracts a pore network using a grain based algorithm |5], which segments the 
pore space into pore bodies and pore throats, each associated with a point x & V, 
a volume, a shape factor G [26] , and an inscribed radius r [32] . One such network is 
visualized in Fig.|3| The distance between the two pore bodies tl and t2 connected 
by a pore throat t3 is divided into parts: parts ki and lt 2 are associated with each 
pore body, tl and t2, respectively; and part Zo associated with the pore throat. 
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Table 1 Micro-CT and model results. 



Porosity 

Ks 

Tortuosity 

Const. 

Char.length 

Permeability 

Rock 

0 

0s 

[itimf] 

Ts 

T 

Cs 

Lhipm] 

Tc [Atm] 

[(At"i) ] 

[mD] 

mCTa 

0.081 

0.048 

0.86 

0.347 

0.373 

107.97 

78.62 

10.4 

0.041 

41.72 

mCTb 

0.128 

0.114 

4.569 

0.41 

0.425 

28.41 

78.5 

15.1 

0.521 

528.06 

mCTc 

0.176 

0.166 

13.165 

0.449 

0.462 

19.91 

102.06 

21.3 

2.192 

2220.86 

mCTd 

0.21 

0.2 

14.078 

0.447 

0.464 

20.49 

107.39 

20.3 

2.811 

2848.57 

a 

0.086 

0.056 

0.766 

0.337 

0.355 

105.3 

75.28 

9.1 

0.043 

43.44 

b 

0.101 

0.079 

1.865 

0.363 

0.381 

58.97 

81.73 

13.2 

0.147 

149.22 

c 

0.125 

0.111 

3.841 

0.387 

0.405 

40.51 

91.08 

15.3 

0.427 

432.61 

d 

0.153 

0.143 

9.573 

0.431 

0.451 

26 

103.6 

18.8 

1.371 

1388.87 

e 

0.176 

0.168 

13.987 

0.442 

0.462 

22.75 

114.03 

22.9 

2.345 

2376.21 

f 

0.206 

0.198 

22.428 

0.461 

0.48 

18.36 

124.6 

26.3 

4.449 

4507.45 

g 

0.245 

0.237 

33.435 

0.467 

0.487 

15.81 

139.26 

30.8 

7.935 

8040.4 


Following Eq. (18) in Ref. [^, the hydraulic conductance of each part is ap¬ 
proximated by 


3r^ 

80pG' 


(45) 


The effective hydraulic conductance between two pore bodies tl and t2 connected 
by a pore throat tS is taken as the length-weighted harmonic average of the three 
parts 


gt = h 


f ill _|_ ill _|_ 

\gti gtz gt 2 J 


(46) 


where It = hi + hz + hz- Assuming a Hagen-Poiseuille type relation between the 
fluid discharge Qti 2 from pore tl to pore t2 and the head gradient ^ = 
{ht 2 — hti)/lt for each pore throat t, we have 


Qtl,2 


= -gt 


pg{ht 2 - hti) 
k 


(47) 


where hu is the piezometric head associated with the pore body ti. 

The network model can now be viewed as a resistor network analog, with a 
one-to-one correspondence between the pore throats in the porous medium and the 
resistors in the resistor network analog, and also there is a one-to-one correspon¬ 
dence between the pore bodies and the network nodes. Each pore throat (resistor) 
t is given a conductance gt/k- Let hi be the piezometric head corresponding to 
pore body (node) i, and the pore throats (resistors) connected to pore 

body i. We then solve for hi such that 


E= f;o,., = o, 


i=i 


(48) 


1=1 


where we have fixed piezometric head at the inlet and outlet boundaries. 
The network volume with non-zero head gradient can be calculated as 


Us - ^ 

hti ^ht2 


y^+y^+vtz 

a?ti at2 


( 49 ) 
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Fig. 4 Plot showing the correspondence between porosity and (ps for the rock models and 
the micro-CT (mCT) data, together with a linear fit to the model data. 


The values for porosity (p and (ps = fPslV for the network representations of our 
micro-CT images and models are reported in Table [1] In Fig. [3] we have plotted 
porosity (p versus (ps. A linear fit to the model data plotted in Fig. |3] gives the 
correspondence 

(Ps{(p) = l.l26(cP- 0.030). (50) 

The fluid velocity inside the network elements is not resolved; we therefore 
threat the fluid velocity as constant inside each network element. The average fluid 
velocity for part tl and t2 are given by uu = gtiPg\hti — hit\ati/Vti where i = 1 , 2 , 
while for part t3 it is uts = gt 3 Pg\hit — h 2 t\/Vt 3 . Here hu and /i 2 t are piezometric 
heads such that gtipg\hti - hit\/lti = gt3pg\hit - h2t\/lt3 = gt2pg\h2t - ht2\/lt2. 
The fluid velocity u is in the opposite direction of the gradient of the head, i.e. 
V/i-u = —||Vh||u. The local permeability factors for the sections as given 

by Eq. (0) , are: 


— 


gti{pg{hti - hit))‘^ati 

hi^ti 



for i = 1, 2, and 


^ _ gt3{pg{hit - h2t))‘^ ( As\^ 

lt3Vt3 \Ah) 


This enables the calculation of the effective permeability factor as the volume 
average of these local contributions: 


f^s 




= 7r E (—-*i + —+ 

j/o \ rv+1 rv+o / 


Vt2 


ottl 


Ott2 


(51) 


The results are reported in Table [T] 

The effective permeability factor Ks versus porosity cp is plotted in Fig. (5] The 
function 


Ks{(p) = 1181 (<(>- 0 . 054 )^■^^ 


( 52 ) 
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Fig. 5 Plot showing effective permeabil¬ 
ity factor Ks versus porosity (p for the 
Fontainebleau rock models and the micro- 
CT (mCT) data, together with a fit to the 
data. 



Fig. 6 Plot showing porosity 0 versus 
permeability k for the Fontainebleau rock 
models and micro-CT (mCT) data, to¬ 
gether with Eq. lisa describing the corre¬ 
lation. 


is included as a fit to the model data. For 0 = 0.054 we have Ats(0) = 0, which is 
interpreted as percolation threshold for this sandstone m- 

The permeability k = Ks4>s, as given by Eq. m, is also listed in Table [1] Cal¬ 
culated porosity and permeability for the network representations of Fontainebleau 
sandstone are plotted in Fig. [6l Combining Eqs. (11011 . (I50II and (15211 . we have: 

fc(^) = Ks{<p)<t>s{(l>) = 1181(0 - 0.054)^-^^l.126(0 - 0.030) 

= 1329(0 - O.O54)^-^^(0 - 0.030). (53) 


This function is also included in Fig. |6l and provides a derived porosity- 
permeability relationship for the Fontainebleau samples. 

We discretized the volume Sis into a disjoint union U5, where 5 is a subvolume 
of a single series of pore throats tg = {t}, with the first throat connected to an 
inlet boundary and the last connected to an outlet boundary. Each 5 transports a 
constant discharge Qg, and y) Qg = Q. The discretization U5 is a simplification of 
the streamlines, similar to the concept of a bundle of capillary tube model. In the 
network representation for the Fontainebleau sandstone the discretization Qs = Ll5 
is dependent on the fluid flow across the network nodes, however different ways 
of tracing streamlines across the network nodes were tested, yielding comparable 
results. 

For the sections ti,t 2 and ts of constant hydraulic conductance as associated 
with pore throat t G tg, we have the corresponding parts {5ti} of 5. Each volume 
Sti C Vti then has the associated length 

The three sections have 


Bti — — 
Bt3 = — 


u 

pgVh 

u 

pgVh 


1 ^2, and 
yti 

Vt3 ■ 


Eor each 5 in 17s = Ll5 we calculated 

B{S) = — StiBti + 5(213(2 + 5(313(3. 
tets 


( 54 ) 
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We separately calculated the constriction factor 


C{S) = 


Ah 


E 


HI 


+ 


H 


+ 


72 

HZ 


\ht2-h2t\ \h^t-h2t\)' 

and the tortnosity 


-(5) = 


As 


E 


(55) 


(56) 


iGts 


We then obtain the characteristic length, constriction factor and tortuosity for 
the volume Qs as: 



(57) 

Cs = ^Y.QsC{S), 

(58) 


(59) 


The calculated values are reported in Table[T] We see that Tg L\/{^Cs) = which 
is consistent with Eq. dSl). 

For comparison we calculated the tortuosity r as 


As 

^ Ylt Qt^t 


(60) 


where Qt is the discharge through pore throat t wm- Note that the values for 
Ts and T are similar, however ts < r. We also calculated the critical pore ra¬ 
dius Cc corresponding to the smallest network element radius of the set of largest 
network elements that percolate through the network m, where the radius of a 
network element is given by r •^3/(10G' tt). The values are reported in Table[T] Such 
characteristic length scales are seen to be signihcantly lower than the hydraulic 
characteristic lengths calculated according to Ref. prz] for the Fontainebleau 
networks, included in Fig. [T] 

In Fig. 0 we have plotted porosity 0 versus the characteristic length squared 
Z/^, the critical pore radius rc, together with a functional ht T^(<(>) = exp(7.650(()-|- 
8.073). For high porosities, i.e. rocks with little cementation and then larger pores, 
we have larger characteristic length Lf^ than for low porosities. The characteris¬ 
tic lengths I/^ for the micro-CT images scatter around the trend given by the 
simulated rock models. 

In Fig. |8] we have plotted porosity (f> versus both the tortuosity Ts and the in¬ 
verse constriction factor CJ. For high porosities we have less permeability reduc¬ 
tion due to both tortuosity and constriction, compared to low porosities. Cemen¬ 
tation increases the ratio between pore body and pore throat cross-sectional area, 
which yields a larger fluctuation in pore size along the streamlines and therefore 
a higher constriction factor. When cementation blocks pore throats completely, it 
increases the length of the streamlines and reduces the tortuosity ts . 

A function r|((^) = 0.108 ln(<()) -1-0.380 gives a visual match to the tortuosity of 
the Fontainebleau sandstone models, while the constriction factor of the models 
follows a trend CT(<(>) = 0.342((7i —0.051). Both functions are also plotted in Fig.P 
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models 

mCT 

exp(7.650*x+8.073) 
Tg models 
r^mCT 


Fig. 7 Plot showing porosity 0 versus characteristic length squared and critical pore radius 
squared for the Fontainebleau rock models and the micro-CT (mCT) data. 



Porosity (fraction) 


o 

§ + models 

o □ mCT 

S X Cs models 

_c * 

O CsmCT 

- 0.108*ln(x)+0.380 

.. 0.342*(x-0.051) 


Fig. 8 Plot showing porosity (j> versus both tortuosity squared and inverse constriction Cg 
for the Fontainebleau rock models and the micro-CT (mCT) data. 


The calculated tortuosity and inverse constriction of the micro-CT images 
follow the trend given by the models. 

The calculated pore structure descriptors for each Fontainebleau sandstone 
sample reveal a strong functional relation with respect to porosity. This is desirable 
for sensitive descriptors. The functional trends display a non-trivial behavior at 
the percolation threshold derived in Eq. ([5^: The tortuosity and characteristic 
length indicate a non-zero value at the percolation threshold, while the inverse 
constriction factor tends to zero at the percolation threshold. 


7 Conclusion 

In this work we have fundamentally described and calculated the permeability in 
porous media. The permeability fc of a porous medium is equal to Ks4>s- Here the 
effective porosity (f)s is the fractional volume conducting flow from inlet to outlet. 
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An effective permeability factor Ks is given by the volnme-weighted average of the 
microscopic permeability factors 

/ \ 2 

K = —apqVh ■ u ( - 7 ^ I . 

\pgAhJ 

This microscopic permeability factor n relates the local contribution of the pore 
structure to effectiveness of the pore space to conduct fluid flow Ks . 

We have shown that Ks = TgBs/Cs = TgL\/[9,Cs), where the effective pore 
radius in the porous medium is described by the characteristic length fluctu¬ 
ation in local hydraulic radii is described by the constriction factor Cg, and the 
effective length of the streamlines is described by the tortuosity Tg. These charac¬ 
teristic length, constriction factor and tortuosity are direction dependent intrinsic 
descriptors of the pore structure. Their directional dependence leads to anisotropy 
of the permeability, i.e., the tensorial form of the permeability. 

We have shown that our methodology reproduces results for Hagen-Poiseuille 
flow in tubes. It is also applied to a natural porous medium given by a pore net¬ 
work representation of Fontainebleau sandstone, where we show how the distinct 
contributions to the permeability from characteristic length, constriction and tor¬ 
tuosity correlate with porosity. As long as the flow and piezometric head field can 
be obtained, this methodology is applicable to any porous medium. 

This work demonstrates how the permeability can be related to porosity, in 
the sense of Kozeny-Carman, through fundamental and measurable descriptors 
of the pore structure. Such derived physical relation between permeability and 
porosity from detailed pore structure information leads to a better fundamental 
understanding of structure-property relations in porous media. 
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